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1 Introduction 

Classically, the properties and dynamics of materials systems can, in principle, be de- 
scribed with the help of their energy landscape. At each given moment in time, we 
associate with the current configuration of N particles a point in a 3 N- dimensional space 
and an energy. The dynamics of the system is reflected in the movement of this point on 
the energy landscape. However, this energy landscape is highly complicated, and trying 
to analytically compute the dynamics is not possible. Thus we are thrown back to the use 
of numerical simulations, which gravely limit the time window over which we are able to 
perform accurate calculations. Furthermore, it can be rather difficult to extract physically 
relevant results from large multi-particle simulations. 

However, the systems often possess features that can aid our quest. Many properties of 
a given system do not depend on time or change according to simple (phenomenological) 
laws, at least within the accuracy of our measurements. This allows us to understand the 
behavior of the system quite well by using simple (phenomenological) models. Frequently, 
the underlying reason for this success is the applicability of separation of time scales 
arguments that one can use to distill the important degrees of freedom relevant on short 
or long time scales. 

Some of the typical processes, each possessing characteristic time scales, relevant in 
materials systems are e.g. 
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• Relaxation of electronic degrees of freedom 

• Lattice vibrations (phonons) 

• Breaking/formation of chemical bonds 

• Changes in the structure of sub-systems, switching between modifications, growth 
of nuclei 

• Grain boundary movement / annealing 

• Diffusional transport, and relaxation to equilibrium due to diffusion 

• Decay of non-equilibrium phases 

• Appearance of fluctuations of certain size (formation of nuclei) 

• Transitions among double (multi) well systems 

• Equilibration of various physical properties 

• Changes in the shapes of sub-systems 

• Interactions among sub-systems 

• . . . 

After a short general discussion concerning simplifications of dynamics from the point 
of view of time scales, we will present examples of relevant time scales and their ap- 
plication to modeling for two systems, amorphous silicon boron nitride, a-Si3B 3 N 7 , and 
homogeneous sinters. 

2 Time scales and modeling 

Materials systems consisting of many atoms, are characterized by a very large number 
of degrees of freedom (xt) = (x± : . . . , xn) (N rs 10 23 ). As a consequence, a complete 
description of their dynamical properties by solving the appropriate quantum mechanical 
equations over the many time scales of interest constitutes an impossible task. 

However, once one is willing to restrict oneself to particular questions and properties 
concerning such a system, it is in many cases feasible to construct a simplified model 
of the system whose dynamics is amenable to numerical or (ideally) analytical solution. 
Technically, such simplifications can be achieved by a number of approaches: 
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• Direct reduction in the number of degrees of freedom: 



(x u ...,x N ) -> (x 1 ,x 2 ). (1) 

Such an approach is particularly useful, if we are interested in a limited time range, 
over which we want to compute the dynamics, which only depends on the variables 
(xi,X2) in that time-range. If (xi,x 2 ) and (x 3 , . . . , xjv) are fast and slowly varying 
variables, respectively, we can replace (x 3 , . . . ,xn) by some constant initial value 
(x 3 , . . . , Xjv) ~ (23(0), . . . , Xat(O)). On the other hand, if (xi, x 2 ) vary slowly, while 
(x 3 , . . . ,xn) fast fluctuate around some average value, it is appropriate to set the 
latter equal to this average (x 3 , . . . , xjy) ~ (#3 , ■ ■ ■ , x n^)- 

Mathematically, this kind of analysis can be performed in a formal way, and we 
refer to the literature on differential equations for more detail 

• Simplification of the dynamics, up to the point where analytically solvable equations 
for the observables of interest, [A s ], result (on a restricted time scale), 

[A s (xt, x N )] -> [A s (xx,..., x N )]. (2) 

This is a very common approach, which consists of applying physically reasonable 
approximations not to the set of variables but to the dynamical equations instead. 

Often, one can combine this approach with the one mentioned above, and, in addi- 
tion, reduce the number of variables one considers, [A a (sci, . . . , x^)] ~ [A s (xi, x 2 )\ 

• A related approach consists of focussing the analysis on a small number of quantities 
of interest [A B (x\, . . . , x^)], which might exhibit a simpler dynamics, possibly de- 
pending on fewer variables, than the full system with all observables . . . , xjy)], 
s = 1, . . . , m: 

[A s (xi, . . .,x N )] -> [Ai(xi, . . .,x N ),A 2 (xi, . . .,x N )) -> [A^xi, x 2 ), A 2 (x!, x 2 )). (3) 

• Finally, and this constitutes perhaps the crowning achievement in a separation of 
time scales approach, one can set up a hierarchy of dynamical systems in a stepping- 
stone procedure. Starting with the full system with all degrees of freedom on the 
fastest time scale (1), one successively reduces the number of degrees of 
freedom until one derives a (phenomenological) set of equations for those observables 
that are relevant for the dynamical description on the longest times, up to infinity. 

One should note that during this procedure, we might employ different sets of ob- 
servables [A^] and variables (x®) for the description of the dynamics on different 
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time scales (i). Furthermore, in many cases, the dynamics on the longer time scale 
will depend on parameters a that incorporate the outcome of the fast dynamics, 
such as a diffusion constant computed from an atomic simulation. In general, these 
parameters a will parametrically depend on the values of the slow variables as pre- 
scribed for the fast-scale calculation: 

[4 i) (4 ) )]-[4 m) (4 m) ;«(4 i+1) ))l- (4) 



In many cases, a rigorous mathematical analysis of the approximations involved is not 
possible or insufficient, and needs to be supplemented by physical reasoning born from ex- 
perience. For materials systems, three general concepts have turned out to be particularly 
helpful, partial (themodynamical) equilibrium, local (thermodynamical) equilibrium, and 
stationary state approximations. 

• Partial equilibrium describes a situation, where the time scale of some degrees of 
freedom is so short that they de-couple from the rest of the system in the following 
sense: Before the slow degrees of freedom change appreciably, the fast ones reach 
thermodynamical (statistical) equilibrium among themselves. This allows two kinds 
of simplifications: 

On very short time scales, we can study the relaxation behavior of the fast degrees 
of freedom while keeping the slow degrees of freedom fixed. 

Conversely, one can study the long-term dynamics of the slow degrees of freedom 
by replacing the fast degrees of freedom by their statistical averages. Such averages 
might depend on the current values of the slow variables, and thus require the 
introduction of additional (phenomenological) terms in the dynamical equations for 
the slow variables on these variables via 

^ = F(xi, x 2 , x 3 , . . . , x N ), i = 1, . . . , N -> 

dx ' 

—j- « F(x u x 2 , x a 3 v (x u x 2 ),..., x%(xi, x 2 )), i = 1, 2. (5) 

Furthermore, in many situations, one needs to include the fluctuations in the fast 
variables (around the respective averages) in the dynamics for the slow variables. 
One common way to achieve this is by adding the fluctuations as a stochastic force 
acting on the slow degrees of freedom. 

• A second, related type of simplifications derives from the concept of local (thermo- 
dynamic) equilibrium. This refers to the equilibration of all (!) degrees of freedom 
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on the (3N-dimensional) energy landscape of the system. As a consequence, the 
energy landscape can be divided into a set of locally ergodic regions 7Z, where the 
system can be considered being in (metastable) thermodynamic equilibrium |3], |J. 
However, in contrast to the partial equilibrium, this state of affairs only persists for 
a limited time, until the escape time r esc from the locally equilibrated regions 7Z 
has been reached. Thus, the local ergodicity of 7Z is only present on a certain time 
scale, 

7~eq ^ ^obs ^eso (6) 

where r eq indicates the time required for the regions 7Z to equilibrate. Mathemat- 
ically, r eq is defined as the smallest time t' for which the difference between time- 
and ensemble average of observables in 71 falls below a certain prescribed accuracy 
a for all t > t', 

\(0) ens -{0) t \<a Vt>t>>r eq . (7) 

For t obs > r esc , the regions 7Z lose their special relevance. Usually, some larger regions 
TV take their place, which are locally ergodic for times between r eq (7Z') > r esc (7Z) 
and T esc (W). 

From the point of view of modeling the dynamics, the relevant configuration space 
can be viewed as being essentially discretized for time scales t b s > T eq (7Z). Here, 
the regions TZ correspond to points on a metagraph (decribed via variables (yj)), on 
which the dynamics of the (possibly new) observables proceeds, e.g.: 

[AV(x 1 ,... 1 x N )]^[A^(y 1 ,y 2 )]. (8) 

In certain situations, one can treat (yj) as if it were a set of continuous variables. 
Describing the dynamics commonly requires developing phenomenological determin- 
istic or stochastic models for (yj) and/or 

In many situations, we encounter a combination of partial and local equilibrium. 
This commonly occurs, if we are dealing with local spatial equilibrium, where the 
system can be modeled as a union of spatially separate pieces, with each being in 
local equilibrium. 

• A special instance of local/partial equilibration that allows the use of separation of 
time scales arguments is the so-called steady-state-situation. Here, external forces 
are applied to the system, such that essentially constant matter/energy flows are 
created in the system. Thus, the system never reaches a "traditional" thermody- 
namic equilibrium, since e.g. the temperature is no longer constant throughout the 
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system. However, we might be able to establish a time-independent temperature 
distribution throughout the material, which would correspond to a locally/partially 
equilibrated state under certain boundary conditions. This allows us to separate 
slow and fast degrees of freedom, where the fast ones equilibrate within a short 
time, and the relevant dynamics of the system on long time scales can be modeled 
using only a small number of variables. 

3 Amorphous Silicon Boron Nitride a-Si3B3N7 

In an amorphous system such as a-Si3B 3 N7, processes take place over a range of very 
different time scales. Many of these are also present in typical crystalline solids, but the 
fact that amorphous systems do not reside in the region around the global (or 'for all 
practical purposes' global) minimum of the energy landscape leads to a steady drift in 
configuration space on all time scales. Although local ergodicity is usually present on 
certain time scales, we always find an escape time beyond which the system must be 
treated as unstable, a process which is commonly denoted as aging. In contrast, systems 
that explore the neighborhood of the global minimum (or any minimum with a similar 
barrier structure representing a stable crystalline modification) of the landscape usually 
possess a limiting time scale, beyond which the system can be treated as globally ergodic, 
i.e. as being in global equilibrium. 

Nevertheless, simplified modeling of amorphous systems using separation of time scales 
is clearly possible. In the following, we are going to list some of the most important time 
scales (with estimates of their value) and discuss their use in modeling the dynamics of 
a-Si3B3N7. 

• Movement and relaxation of electronic degrees of freedom (r rs 10 -15 sec). 

This is the time scale on which electronic, magnetic or optical properties, need to be 
addressed. From the point of view of a stepping stone approach for the longer time 
dynamics, one would attempt to employ ab initio calculations to derive empirical 
potentials that can be used in subsequent steps when performing Monte Carlo or 
Molecular Dynamics calculations. For a-SiaBaNy, several such potentials have been 
constructed m, which have been employed in the calculations discussed below. 

• Lattice vibrations (10 -9 sec > r mb > 10~ 14 sec) 

We note that if we perform simulations of the system in the picosecond region, 
we are able to study the vibrations of the network, because the time scales on 
which the network itself changes (t bsp , see below) are longer by far, at least at low 
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Figure 1: Typical vibrational spectrum along trajectories for five different temperatures 
between 25 K and 4000 K. 'Negative' frequencies u correspond to imaginary modes. 



temperatures. This can be seen in the vibrational spectra along the trajectories (c.f. 
fig. [l]). There, we observe that for low temperatures, only a very small percentage 
indicate imaginary modes, i.e., that the trajectory is in a neighborhood of a saddle 
point in the corresponding eigendirection. Clearly, most of the time we are in the 
process of vibrating around a local minimum. Only at high temperatures close to 
or above the glass transition point (T c , see below) do we encounter a substantial 
fraction of imaginary modes. 

Encouraged by the separation of time scales argument, one can study the spectrum 
of real eigenvalues of the minima located along the actual trajectories (c.f. fig. 0). 
Here, we find that the energy landscape of the liquid (T > T c ) and the glassy 
(T < T c ) state are very similar with regard to the local minima lying below. The 
general shape of the frequency spectrum of all minima is nearly identical, and one 
only observes a slight change towards a higher average frequency at minima with 
lower potential energy. Since the deeper-lying local minima are encountered more 
frequently in the glassy state, there should therefore occur a slight softening of 
the vibrations when passing through the glass transition region upon increasing 
the temperature. But we note that this effect is not large enough to account for 
the peak in the specific heat (c.f. fig. |3|) observed in this temperature range 0. 
Regarding this figure, one should note that the contribution of the kinetic energy 
(3/2/cb) has already been subtracted. For T — > oo, Cy drops below 3/2fcg because 
the anharmonic part of the potential becomes dominant and the potential goes to 
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Figure 2: Typical vibrational spectrum for local minima of the energy landscape observed 
along trajectories for five different temperatures between 25 K and 4000 K. 

3.5 
3 

_ 2.5 

> 
O 

2 
1.5 
1 

2000 4000 6000 8000 10000 

T[K] 

Figure 3: Specific heat per atom calculated via Cy = A(E)/AT. 

zero for long distances. For T < T c , Cy falls below 3/2kB because of aging effects 
involved in the numerical computation of (E)(T) for neighboring temperatures T 12 
(AT = T 2 — Ti). This effect is a consequence of the different rates of decrease m in 
minimum energies with temperature shown in fig. |5], which have a maximum around 
1750 K. 

• Time scales associated with the transport of energy 

A typical time scale of energy transport that is relevant during molecular dynamics 
simulations is given by: t$ oc L/vg ~ 10~ n — 10 -13 sec, where L is the length 
of the simulation cell (typically, L « 10 — 100 A), and v$ is the speed of sound 
(typically, v$ ~ 10 3 — 10 4 m/sec). Turning this relation around, we can estimate 
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the speed of sound in our system by studying the time evolution of a heat pulse 
that is excited along one layer of the system at a given time. From the speed 
with which the energy spike travels through the system, we can deduce a value 
for vs ~ 3 — 4 x 10 3 m/sec for the energy nowH, a value which lies in the range 
of the speeds of sound found experimentally in crystalline SiaN 4 (v$ »s 1.1 x 10 4 
m/sec) and BN (v$ ~ 7 x 10 3 m/sec). Simulations of energy transport in SiaN 4 and 
BN, using the same interaction potential, yield velocities vs ~ 7.5 x 10 3 m/sec and 
vs ~ 8 x 10 3 m/sec, respectively ||. 

This transport time Ts competes with the time needed to dissipate the energy 
throughout the system, mostly due to phonon-phonon-interactions. This dissipation 
time also can be estimated from the heat pulse experiment || , to yield Tdi SS ~ 50 fs. 
One should note that this value strongly increases with the degree of crystallinity in 
the system (similar observations were made for amorphous Se when compared with 
crystalline Se|J). 

One important application of this analysis is that we can simulate and measure the 
heat conductivity of a-Si3B 3 N 7 by generating a steady-state heat flow through the 
system. To achieve this, we perform the following simulation procedure 0: 

1. Equilibrate the system at constant temperature for a sufficiently long time T wa u 
(see below) to keep aging to a minimum for the remainder of the experiment. 

2. Increase the temperature along a layer of atoms located at L/4 in the simulation 
cell by AT = 0.1 x T, and decrease it be the same amount along a layer at 
3L/4. 

3. One now lets the system relax on the time scale of the dissipation time. By 
then a temperature gradient is established, and a constant flow of heat into 
and out of the system is present. 

4. We maintain this steady-state situation for an observation time t Q b s < T wa u, 
and measure the average heat flow from the hot to the cold plate. 

5. This information can be used as input into the standard phenomenological heat 
transport law, 

K = JjL = m.±!l (9) 

|VT| AAt 2 AT' v ; 

to deduce the heat conductivity k, where \j\ = (AE)/(AAt) is the average 
heat flow through a cross section A, and |VT| = 4AT/L is the temperature 
gradient. 
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Figure 4: Heat conductivity k as function of temperature T. 

6. In a final step, one repeats this procedure for larger simulation cells, in order 
to take the limit of k for infinite cell size (L — > oo). 

Figure £| shows the extrapolated heat conductivity as a function of temperature. We 
note that this value, k « 4 W/mK, shows only a small variation with temperature ||, 
similar to the results obtained for amorphous Se.||] Since no experimental value 
has been measured up to now for a-Si3B 3 N 7 , we estimate the reliability of the 
procedure by studying the corresponding crystalline binary systems. One finds that 
the experimental value exceeds the one from the simulations by about a factor 3 — 5, 
similar to what we had found in the case of Selen||. Thus, one would conclude that 
the true value of k should be about n exp « 15 W/mK. 

• Aging, Dissipation, and Response 

Unless interested in relaxation processes as such, one usually waits until the system 
being simulated has relaxed into a state of partial or local equilibrium. Several 
time scales are of particular interest in this context: The time(s) Ti n u necessary 
to "forget" various aspects of the initial conditions, and the time(s) r resp needed 
to react to the application of an external force. In the latter case, we might be 
dealing with either constant or time-varying, in particular sinusoidal, forces. We 
note that r resp usually implies a relatively weak external disturbance of a system 
that is already equilibrated. In contrast, T init refers to a relaxation from a possibly 
highly non-equilibrium initial condition of our simulations. 
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Up to this point, these considerations are applicable to both crystalline and amor- 
phous solids. But when simulating amorphous systems, one has to take into account 
that the region of the energy landscape we explore is not centered around a (global) 
minimum that is separated from similarly deep minima by insurmountable (on the 
time scale of the simulation) barriers. Instead, a permanent drift towards lower and 
lower energy minima takes place. Although the barriers needed to be crossed to 
reach the lower-lying minima constantly increase, this only leads to a logarithmic 
decrease in the speed with which deeper regions become accessible but does not stop 
this drift. 

As a consequence, we are dealing with a so-called aging process: For observation 
times short compared to the time T wa u the system has been allowed to relax freely, 
the system behaves as if it were in thermal equilibrium, while for times exceeding 
T wait non-equilibrium behavior sets in. This can be seen by studying various two- 
time correlation functions, such as the two-time energy-energy average 

J.U 4- T"i\ \E{t W ) ■ E(t w + t b s )) ens M v 

<f)(t w , t obs ; T) = (10) 

or the closely related two-time autocorrelation function 

Ce^u], tobs] T) = (E(t w ) ■ E(t w + t b s )) ens — (E(t w )) ens ■ (E(t w + t bs)) ens, (H) 

which in a (quasi)-equilibrium regime equals a generalized standard equilibrium 
specific heat kBT 2 Cv(t w ,t b s ',T). 

In figure |5], we show the potential energy E po t of local minima along MC-trajectories 
as function of time \n.t b s for several temperatures averaged over many trajectories ||. 



One clearly notices that after Tj n # ~ 10 3 MCC (Monte-Carlo cycles) the behavior 
changes and a clear downward drift of the energy occurs which is proportional to 
Int. This value of T init is also observed in other quantities such as the diffusion 
of atoms. For t > r init , the 'ballistic' phase appears to be over, and vibrations 
about local minima together with slow hopping from minimum to minimum on the 
landscape best describe the dynamics in configuration space. 

This drift on the landscape leads to aging as can be seen in figure where we show 
the development of the two-time energy-energy average as a function of waiting 
time|7[]. The departure from the (quasi-) equilibrium value (= 1) for t obs > r wait is 
clearly observable. 
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Figure 5: Average potential energy E pot of local minima along trajectories as function of 
time lnt b s for several temperatures. The inset shows the negative slopes m as function 
of temperature. 
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waiting times at T = 1250 K for an ensemble of 100 runs. 
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Figure 7: Bond survival probability BSP as function of time for several temperatures. 

• Change of the Network 

The major events that take place during the aging process are associated with 
changes in the topology of the network. One quantity that reflects this behavior is 
the bond survival probability BSP, which is defined as the probability as a function 
of time that a bond that was present at time t\ is still present at time t<i (without 
having been broken in-between). Figure [?| for the Si-N bonds shows that for high 
temperatures, where the system is essentially fluid, BSP(T ini t,t) follows a stretched 
exponential curve, BSP(t) oc exp(— ) ); the curves for B-N bonds exhibit 
the same functional behavior. A similar behavior is found for low temperatures, 
although here one cannot observe the decrease to zero - on the time scale of the 
simulations, most of the bonds remained unbroken. 

The characteristic time scale is tbsp, which increases exponentially once we go 
below a critical temperature T c m 2000 — 2500 K. Quite generally, the behavior of 
a-SisBsNy is comparable to that of a glass when crossing this temperature region: 
A freezing-in of the structure, associated with a peak in the specific heat and an 
exponential increase of the viscosity 0. For temperatures below T c , we observe some 
weak aging in the bond survival probability, i.e. the network becomes more resistant 
to change as time passes. This is not surprising considering the downward drift in 
the potential energy. 

Nevertheless, the existence of a typical time scale for bond-breaking allows us to 
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apply separation of time scales arguments to the modeling of the system. For 
tobs > tbsp ( t bsp > 10~ 7 sec for T < 500 K), we can dispense of the vibrational 
degrees of freedom, and choose the topology of the system as the relevant quantity. 
It is then possible, to replace the continuous configuration space by a metagraph, 
where each node corresponds to a (possibly stretched) network topology, and the 
edges of this graph correspond to the possible bond-breaking and bond-stretching 
moves of the atoms. In order to assign realistic energies to these configurations, one 
would place the network on a lattice with an appropriately chosen lattice-spacing, 
and study the resulting energy landscape, and model the dynamics as a stochastic 
process on the metagraph. 

Such work has been performed in the past for a number of discrete complex systems, 
such as the TSP-landscapeJIIJ, chain-polymers [pi]] or one-component 2d- and 3d- 
lattice networks [[12], [13|]. While such calculations have not yet been attempted for 
a-Si3B3N7, it would be a natural further step to take, in order to model dynamical 
properties of this system for very long times. 



4 Homogeneous Sintering 

Once we leave the realm of physical phenomena that can be described on picosecond time 
scales, we are even more urgently in need of physically reasonable simplifications, if we 
want to model and understand the various processes taking place during the generation of 
a macroscopic material, and its properties on macroscopic time scales. As an example we 
are going to analyze homogeneous sintering, and show how one can use a stepping stone 
approach to model the complete sintering process starting with clusters and ending with 



a dense macroscopic sinter. W% 15 



We begin with a short summary of the sintering process, where for simplicity it is 
assumed that only one type of (metal) atoms is involved. Typically, at the outset a dense 
packing of small particles or atom clusters is prepared, which in many cases has already 
been arranged in a particular macroscopic shape needed as a part in a work-piece. This 
conglomerate is now heated to about 2/3 of the melting temperature of the bulk material, 
leading to an enhanced diffusional motion of the atoms. The total amount of 'liquid' (= 
mobile atoms) present is relatively small, and the necks between the clusters are wetted. 
During this time, new 'liquid' will be produced throughout the sinter at places favoring 
high atomic mobility (e.g. naked surface of small spheres), while some of the 'liquid' will 
solidify in thermodynamic sinks of the material (e.g. in the interior of the necks connecting 
the spheres). At the beginning, the communication among the various 'liquid-like' wetted 
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regions is rather slow, and the shapes of the neck-configurations correspond to locally 
equilibrated arrangements of the mobile atoms. 

The surface tension generated by the wetting process leads to a net-movement of 
the clusters towards each other, in order to minimize the total surface energy, and a 
redistribution of the atoms at the necks and on the surface of the clusters. This process 
competes with the exchange of 'liquid' material among the various pores and necks, which 
would lead to a global equilibration among the wetted regions. At this stage, an important 
role can be played by the heat-flow through the system. This flow is rather complicated 
due to the many melting and re-freezing processes that can, in principle, take place within 
the sinter. 

Once a certain densification of the solid-like skeleton has been achieved, the open 
channels providing for easy surface diffusion-driven communication among the 'liquid-like' 
regions are closed off, and atom transport is now only possible via grain boundary diffusion 
or bulk diffusion. At this stage of the process, it is energetically favorable for the isolated 
pores to become more spherical. Furthermore the annealing of the grain boundaries that 
had been generated at the necks between the spheres and along the closed off channels 
begins. In addition, the smaller pores start drifting, and mergers among them can occur 
at this stage of the process. 

These processes continue into the final stage of sintering, where a coarsening process 
among the surviving (spherical) pores takes place. Here, atoms are exchanged via bulk 
and grain boundary diffusion, leading to a net growth of larger pores at the expense of 
smaller ones. In the ideal case, no pores will be left after a very long time, all "vacancies" 
having diffused to the surface. During the whole process outlined, the sinter is likely to 
have shrunk in volume by 20 - 40 %. Controlling the local macroscopic shrinkage rates is 
therefore an important issue in the efficient production of pre-shaped work-pieces using 
the sintering process. 

When trying to model this complicated process, preferably starting at the atomic scale, 
separation of time scales arguments turn out to be very helpful. Some of the relevant time 
scales for sintering are: 

• Relaxation of electronic degrees of freedom (r el < 10~ 15 sec) 

• Lattice vibrations (10~ 9 > r vib > 10 -14 sec) 

• Achieving local (spatial) equilibrium during wetting (10~ 4 > r eq > 10~ 8 sec) 

• Production of the mobile pool (10~ 2 > r pr > 10~ 6 sec) 
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• Atomic transport (communication time) among locally equilibrated wetted regions 
(1(T 2 > r com > KT 6 sec) 

• Equilibration among wetted regions (10~ 2 > r egens > 1CT 6 sec) 

• Movement of solid particles across several atom-atom distances (1CT 2 > r core > 1CT 6 
sec) 

• Annealing (movement) of grain boundaries (10 5 > r gh > 1CT 1 sec) 

• Rounding of pores (10 3 > r ro > 1CT 2 sec) 

• Equilibration of the atoms within grain boundaries (10 3 > r eqgh > 1(T 4 sec) 

• Time needed for atoms to diffuse between isolated pores (10 5 > r d% ^ > 1CT 3 sec) 

Of course, all these times are strongly dependent on pressure, temperature, connectivity 
of the sinter, type of atoms, etc. 

These time scales suggest that modeling of sintering involves two different classes of 
models, the analysis of the atomic movement in an effective potential, and a continuum 
thermodynamic model that employs the parameters deduced from the atomic simulations 
as input. The first class includes the computation of parameters for an effective potential 
of atom movements, followed by molecular dynamics simulations of the behavior of clusters 
and their fusion. Finally it involves the calculation or estimation of various parameters 
and time scales needed in the thermodynamic model, such as diffusion constants, the time 
scales for movements of grain boundaries, the speed with which pores are closed off and 
become spherical, and the time scale for isolated pores to diffuse through the sinter and 
merge with other pores. 

Here, we are going to focus on the thermodynamic continuum model, where we use 
the separation of time scales to identify the following processes: 

• Fastest process: local wetting of pores and necks. 

• Slower competing processes: melting / re-freezing / transport of mobile material, 
equilibration among the wetted regions, movement of solid nano-particles constitut- 
ing the solid skeleton. 

• Slowest processes: annealing of grain boundaries, drift of isolated pores, and coars- 
ening of isolated pores via emission and absorption of atoms / vacancies. 
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This general concept was implemented with several groups participating. J. Bernholc 
and co-workers at NCSU derived EAM-potentials []IE | , and H.P. Cheng and R.S. Berry 
in Chicago studied the melting and fusion of clusters [IT], TE\. Summarizing their most 
important results, they observed surface melting for clusters and found that a neck grew 
between two clusters via wetting and that a grain boundary was created along the neck, 
which only slowly could be removed. 

The process of sintering several clusters, through the closure of pores, their movement 
and merger, and coarsening via vacancy emission was studied using a very simplified 
lattice-gas type model, in order to be able to bridge the many time scales involved [|r9|, |20[| . 
It was found that the closure of pores did occur relatively fast |19| , and that the shape of 
the clusters and pores exhibited strong fluctuations about the spherical " ideal" minimum 
energy shape. For small pores, merging was the dominant coarsening process for high 
pore density (the density of the pores was essentially given by the overall density of the 
sinter). In contrast, for larger pores, standard Ostwald-ripening occurred except for very 
closely neighboring pores, where large shape fluctations (e.g. a massive elongation of the 
pore) could lead to a direct contact among these pores resulting in a fast merger and 
rounding of the merged pore. |20| 

Besides the separation of time scales, the guiding principle of the thermodynamic 
model is the minimization of the (free) energy. This energy minimization corresponds to 
the determination of the average values of the fast degrees of freedom in partial and/or 
local thermodynamic equilibrium. Regarding the initial configuration of our model of the 
sintering process, we note that the starting sinter is created by densely packing small 
particles, perhaps under application of pressure. Thus, a minimum energy configuration 
corresponds to an as dense as possible packing of these particles. Finding the densest 
packing for e.g. spheres of different sizes is a global optimization problem. But even de- 
termining local minima, i.e. arrangements that can no longer be compressed any further 
by infinitesimal shifts of the spheres and are stable against small to moderately large rear- 
rangements of the spheres, is a non-trivial task. Such arrangements of spheres constitute 



random close packings of spheres, which we have generated by appropriate algorithms. |21 



The next step in our approach requires the wetting of pores and necks. At the be- 
ginning, the optimal location for the mobile pool of atoms is within the necks between 
pairs of fixed unchangeable nearest-neighbor spheres. The minimum energy configuration 
is given by minimizing the surface energy of the two-sphere system for a given amount 
of liquid (Plateau's problem). While this task can be solved analytically in terms of 
elliptic functions, we have found it more useful to directly solve the corresponding Euler- 
Lagrange equation with the appropriate boundary conditions, and tabulate the shapes of 
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the optimal solid+liquid configurations as a function of radii of the spheres, center-center 
distance of the spheres and amount of liquid available. For single necks, the minimum 
energy surface is axisymmetric, exhibits constant mean curvature, and the cross section 
is given by a Delaunay curve. |f2"2", 



The transport of liquid material between the wetted regions can be modelled via net- 
work thermodynamics. [EM], |14|| Here, one describes the net-diffusive transport Jy between 
many sources and sinks i and j containing iVj and Nj atoms, respectively, in analogy to 
an electric network consisting of resistors and capacitors. In the model, the conductivity 
and the voltage differences correspond to the effective diffusivity er^- along the channel or 
grain boundary and the difference in chemical potentials fij(Nj) — fii(Ni) (surface energies 
of the wetted regions), respectively: 

^r = £ J « = E ^At'A^) - m(Ni)) (12) 

We note that this model also applies to the description of the interaction via grain bound- 
ary diffusion of a dense distribution of isolated pores at the beginning of the coarsening 
process. 

The movement of the individual solid spheres follows simple Newton's equations, where 
three types of forces contribute: 

• Force, due to minimization of local surface tension: 

Fs = -a^— BK -h i , (13) 

i=i 1 

where hi is the unit vector along the center-center axis between neighbor spheres 
along neck i, a the surface tension, and Aj(|/tj|) the surface area of the neck, respec- 
tively. 

• Damping force due to the viscosity caused by refreezing of the necks: 

F vis = -T)Vi, (14) 

where r\ is the viscosity. 

• Forces due to steric constraints. 

Depending on temperature, viscous forces or steric constraints tend to dominate, in the 
sense that we can assume that the particles either move with a nearly constant ('final') 
velocity (dominated by viscous forces) or move at each step the maximum distance allowed 



until solid-solid contact is established. [Eq, 14 
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If the changes in the skeleton are slow compared to the transport of material among the 
liquid pools and the production of new liquid material via surface melting of the spheres, 
the network thermodynamics model applies. In the extreme case of very slow addition 
of new material to the liquid pool (low rate of surface melting) and nearly unchanged 
solid skeleton, then a global equilibration among the wetted regions can take place. This 



equilibration constitutes again a global optimization problem|26|], and can be solved e.g. 
using dynamical programming f25| . 

On the other hand, if the liquid production is very slow, but liquid transport is fast 
and the solid particles easily movable, we are dealing with a generalization of Plateau's 
problem, i.e., we search for the minimum surface of a combination of many solid movable 
particles plus a freely deformable liquid pool[[27j . 



The most cumbersome task we are faced with is the determination of the detailed 
shape of the wetted regions once several necks coalesce, but the channels between the 
pores have not yet closed. No analytical solutions are available (the solutions are no 
longer axisymmetric surfaces), and the number of possible multi-sphere+liquid configura- 
tions that need to be considered becomes very large. Therefore, tabulating the numerical 
solutions of the Euler-Lagrange equations does no longer appear to be an efficient ap- 



proach, and one has to use e.g. a numerical minimum surface solver [p8| at this stage 
after every liquid-redistribution step. However, our simple model simulations mentioned 
above indicate that this stage does not last very long, and that the pores rather quickly 
become isolated. Once this has happened, rounding of the pores begins to take place, 
together with a slow drift, both of which are amenable to analytical modeling. Depend- 
ing on the size and the density of the pores, the interaction at this stage will be best 
modeled by network thermodynamics or by models of merging pores analogous to the von 
Smoluchowski-approach. [P0f | 

Once the pores can be considered as essentially immovable, standard coarsening pro- 
cesses will take over. Here, one needs to address the question of bulk vs. grain boundary 
diffusion, and the issue of shielding of the pores by neighboring pores. Clearly, the models 
depend on whether the interaction among the pores proceeds via three- or two-dimensional 
diffusion. Bulk diffusion clearly goes via three dimensions, while the case of grain bound- 
ary diffusion is not as clear-cut. We might regain three-dimensional diffusion, if the grain 
boundaries intersect very frequently in 'zero volume' pores or lines. In such a case, the 
vacancies / atoms within the two-dimensional grain boundary network 'notice' the three- 
dimensionality of the embedding space, and the net-diffusion throughout the solid can 
be treated as if it were three-dimensional, although with a different effective diffusion 
coefficient. 
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If the grain boundary network is rather sparse with each grain boundary connecting 
many pores and most pores intersecting one grain boundary only, however, we can replace 
the three-dimensional sinter by a single two-dimensional plane that intersects all the 
three-dimensional pores. [fH|, [2"9| In this situation, we can solve the equations involved 
(2d-diffusion-, 3d-Fokker-Planck-equation, 3d-conservation of vacancies) employing some 
generalization of the standard procedures. 

Again, separation of time scales comes to our assistance by allowing us to separately 
analyze the process of rounding of the pores to spherical cavities, the establishment of a 
diffusional flow of one-atom-size vacancies between the individual pores and the environ- 
ment at infinity, and also the net exchange of material among pores and the growth of 
the average size of the surviving pores. (Note that at the final stage the dynamics takes 
place on a metagraph, whose nodes correspond to the feasible distributions of clusters of 
various sizes, and whose edges correspond to allowed changes of the cluster distribution.) 

The result is similar to the one derived for bulk-diffusion ||, i.e. the scaled radius x is 
related to the scaled time t' via x oc (t') 1 ^ 3 , but with correction terms: 



(lit -Poj) [l - \{m 2pi W - Po ) ( - 1/3) + 0(f 2 (t> -p )- 2/3 ) 



X 



(15) 



where po and /2 are combinations of several constants characteristic of the sintering ma- 
terial and its state at the beginning of the coarsening process. 

We note that in the standard approach, the pores only interact via a sea of vacancies 
(at infinity), either in the bulk (3d) or in the infinite grain boundary (2d). If the density 
of pores is still quite high, however, this can no longer be assumed to hold, and shielding 
of the pores by neighboring pores becomes important. Again, this problem can be treated 
analytically on various time scales and we find that e.g. the mean cluster size M in d 
dimensions grows as 

M(t) oct d/(4+c °. (16) 

At first sight, this example of a homogeneous sinter might appear to be far removed 
from the field of non-crystalline and partially crystalline materials, which has been the 
topic of this workshop. But from the point of view of the modeler both the general 
techniques employed and even many of the particular algorithms presented can easily be 
adapted to the analysis and modeling of macroscopic dynamical processes in amorphous 
systems, such as the growth of crystals from the amorphous phase or the generation of 
nanocrystalline material with a specific texture. 
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